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Introduction 

The  first  full  year  of  research  for  the  project  entitled  “Bayesian  Hierarchical  Models  (BHM)  to 
Augment  the  Mediterranean  Ocean  Forecast  System  (MFS)”  completed  at  the  end  of  May  2006. 
Project  achievements  have  met  or  exceeded  plans  put  forth  in  the  proposal.  Prof.  Nadia  Pinardi 
(Univ.  Bologna,  INGV)  and  Dr.  Ralph  F.  Milliff  (NWRA/CoRA)  presented  early  results  to  Phys¬ 
ical  Oceanography  Program  Managers  in  a  seminar  at  ONR  Headquarters,  in  Arlington,  in  early 
May.1  This  annual  report  reviews  highlights  from  that  presentation. 


MFS  -  Wind-BHM 


Initial  project  activities  focused  on  the  development  of  a  BHM  for  the  surface  wind  process,  based 
on  the  stochastic  geostrophy  model  first  proposed  by  Royle  et  al.  (1998).  The  MFS-Wind-BHM 
implementation  was  led  by  Prof.  Chris  Wikle  (U.  Missouri,  project  Co-PI).  Several  versions  of  the 
surface  wind  BHM  have  been  implemented  and  tested.  Version  4  has  been  implemented  in  the  first 
MFS  ensemble  data  assimilation  and  forecast  experiments  reviewed  here. 

In  probability  model  notation,  MFS-Wind-BHM  is  given  at  the  lowest  (i.e.  most  general)  level 
in  two  parts;  a  data  stage  distribution  and  a  process  model  stage  distribution.  These  are: 

T  T  T+L 

H[St\Wt,ds]xl[[Ar\Wt,dAWM\Pt,eAP]x  H  [F?\Wu6FW}[F?\Pt,6FP\-  (1) 

t=  1  t= 1  t=T+ 1 


N  x  X\[wt\puew][PtKep\[bt\bt^ob}  (2) 

t.= i 

where  probability  distributions  are  noted  by  square  brackets,  and,  for  example,  the  conditional  dis¬ 
tribution  for  the  wind  process  Wl  given  the  SLP  process  Pl  is  denoted  [W^P*].  In  (1)  the  data 
stage  distributions  for  the  data  assimilation  period  t  —  1,T  include  scatterometer  winds  S,  and 

'CDs  containing  the  powerpoint  slides  used  in  the  seminar  have  been  sent  to  Drs.  Scott  Harper,  Manuel  Fiadeiro 
and  Theresa  Paluszkiewicz,  Physical  Oceanography  Program,  ONR. 
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analysis  surface  wind  Aw  and  SLP  Ap.  During  the  forecast  period  t  =  T,T  +  L,  the  ECMWF 
forecast  winds  ( Fw )  and  SLP  (Fp)  degrade  in  skill  with  increasing  forecast  lead  time.  Following 
Royle  et  al.  (1998),  the  process  model  part  of  (2)  is  based,  in  part,  on  a  stochastic  geostrophy 
model,  wherein  the  surface  wind  distribution  [IF*]  is  conditioned  on  gradients  of  the  SLP  process 
P*  (both  cross-  and  down-gradient  terms).  SLP  fields  are  constructed  by  time-dependent  weight¬ 
ings  bt  of  spatial  structure  functions  (i.e.  eigenvectors).  Time  dependence  in  the  prior  is  prescribed 
in  an  autoregressive  model  at  the  next  higher  level  of  the  BHM  hierarchy  (not  shown). 

Figure  1  depicts  a  snapshot  of  10  realizations  from  the  posterior  distribution  for  the  surface 
wind  process  as  generated  in  the  version  4  implementation  of  MFS-Wind-BHM.  A  cluster  of  10 
blue  surface  wind  vectors  is  depicted  at  every  other  BHM  prediction  grid  location  (prediction  grid 
resolution  is  0.25°)  in  the  western  basin  of  the  Mediterranean  Sea.  A  green  vector  within  each  clus¬ 
ter  depicts  the  ensemble  mean  surface  wind  for  each  location  (i.e.  the  average  speed  and  direction 
of  the  10  realizations).  A  red  vector  (sometimes  obscured)  depicts  the  ECMWF  analysis  vector 
for  the  same  output  time.  When  MFS-Wind-BHM  winds  are  organized  by  a  prevailing  weather 
pattern  (e.g.  the  Sirocco  event  depicted  in  Fig.  1),  the  distributions  of  wind  speeds  and  directions 
are  relatively  tightly  spread  about  the  ensemble  mean  vector.  Conversely,  in  low  wind  speed  events 
(e.g.  see  the  Tyrrhenian  Sea),  the  wind  direction  is  arbitrary  as  reflected  in  the  distributions  that 
span  a  wide  range  of  directions,  all  at  very  low  wind  speeds. 

An  animation  has  been  created,  in  the  format  of  Fig.  1,  for  every  MFS-Wind-BHM  output 
time  (4-times  daily)  spanning  the  data  assimilation  and  forecast  periods  (i.e.  t  =  1  ,T  +  L). 
During  the  assimilation  period,  there  are  instances  when  the  ECMWF  analysis  vector  departs  from 
the  MFS-Wind-BHM  cluster.  This  is  due  to  the  influence  of  S  in  the  data  stage.  Evidence  has 
shown  that  the  scatterometer  winds  do  not  influence  the  ECMWF  analyses  at  the  resolution  of 
the  oceanic  mesoscale  in  the  Mediterranean  Sea  (Milliff  and  Morzel,  2003).  Conversely,  during 
the  forecast  period,  departures  of  the  MFS-Wind-BHM  clusters  from  the  (retrospective)  ECMWF 
analysis  vectors  are  probably  due  to  diminishing  skill  in  the  MFS-Wind-BHM  forecast  (i.e.  for 
lead  times  greater  than  about  3  days). 

The  MFS-Wind-BHM  animation  was  shown  and  discussed  in  the  seminar  at  ONR  HQ. 

MFS  Ensemble  Initial  Condition  and  Forecast  Responses 

The  MFS  data  assimilation  and  forecast  cycle  are  depicted  schematically  in  Figure  2.  Ocean  data 
from  the  multi-platform  MFS  observing  system  are  assimilated  over  a  14  day  period  (i.e.  T  =  14 d). 
Data  sources  include  XBT  profiles  from  Volunteer  Observing  Ships,  CTD  profiles  from  occasional 
directed  research  cruises  (including  a  glider  experiment  in  the  Ionian  Sea),  profile  and  current 
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information  from  ARGO  drifting  buoys,  ocean  and  meteorological  data  from  fixed  moorings  (i.e. 
M3A  buoys),  and  sea-surface  height  (SSH)  anomaly  fields  from  accumulated  altimeter  data.  The 
ocean  data  are  assimilated  in  a  Reduced  Order  Optimal  Interpolation  (ROOI)  procedure  that  uses 
ECMWF  analysis  atmospheric  forcing  data  to  force  tangent  linear  dynamics.  A  10-day  ( L  =  10  d) 
ocean  forecast  is  launched  from  the  end  of  the  assimilation  period.  Surface  forcing  for  the  MFS 
forecast  period  is  taken  from  ECMWF  forecasts. 

For  thirteen  days  after  the  ocean  assimilation  step,  an  ocean  simulation  is  run  to  span  the 
period  from  the  last  ocean  analysis  to  the  present  day  (i.e.  simulations  of  1,  2,  3,  ...  13  days  are 
run  on  second,  third,  fourth,  ...,  twelfth  days  after  the  assimilation).  These  simulations  are  forced 
by  the  new  ECMWF  analysis  data  that  become  available  after  the  initial  forecast  run.  New  10-day 
ocean  forecasts  are  initialized  at  the  end  of  each  simulation.  In  this  way,  10-day  ocean  forecasts  are 
generated  every  day  as  the  standard  MFS  output.  Every  fourteenth  day,  a  new  data  assimilation  step 
is  performed  and  the  MFS  forecast  stems  directly  from  the  initial  conditions,  without  a  simulation 
step.  All  other  days,  the  10-day  ocean  forecast  launches  from  the  assimilation  period  plus  the 
simulations,  where  the  simulation  period  increases  by  1  day,  everyday. 

MFS  Data  Assimilation  Step  with  MFS-Wind-BHM 

In  our  reforecast  experiments,  we  use  10  realizations  from  the  posterior  distribution  for  MFS- 
Wind-BHM  to  drive  both  the  assimilation  step  (t,  =  1,  T )  and  the  10-day  forecast  (t  =  T,T  +  L). 
The  basin  scale  statistics  of  the  10  realizations  are  similar  from  realization  to  realization.  Basin 
average  RMS  difference  with  respect  to  the  ensemble  mean  are  around  0.75  ms-1;  dropping  to 
0.55  ms~l  when  scatterometer  data  are  available  during  the  assimilation  step. 

Importantly,  the  ocean  response  to  the  MFS-Wind-BHM  realizations  during  the  assimilation 
step  concentrates  the  largest  RMS  differences  in  mesoscale  eddy  structures.  Figure  3a  is  a  map  of 
the  initial  condition  (i.e.  t  =  lAd)  spread  in  SSH.  Figure  3b  is  the  initial  condition  spread  in  sea- 
surface  temperature  (SST).  The  mesoscale  eddies  are  the  dominant  features  of  the  initial  condition 
spread. 

This  is  an  ideal  response  for  the  MFS  data  assimilation  system,  and  it  is  a  direct  result  of  having 
forced  the  assimilation  step  with  surface  wind  realizations  from  the  MFS-Wind-BHM.  Initial  con¬ 
dition  spread  is  a  measure  of  the  uncertainty  in  the  assimilation  procedure;  uncertainty  that  stems 
from  both  observational  and  model  errors.  Figure  1  depicts  uncertainty  in  the  surface  forcing  that 
is  sensitive  to  the  region  and  meteorological  scenario  captured  in  the  observations  and  analyses. 
Figure  3  implies  that  subtle  changes  in  speed  and  direction  of  the  surface  forcing  drive  subtle  dif¬ 
ferences  in  Ekman  pumping  and  suction  at  the  surface.  These  differences  accumulate  in  mesoscale 
eddies  that  are  the  most  uncertain  features  of  the  Mediterranean  circulation.  An  initial  condition 
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perturbation  method  that  concentrates  uncertainty  in  the  mesoscale  eddies  is  therefore  optimal. 

MFS  Forecast  Step  with  MFS-Wind-BHM 

A  major  goal  of  ensemble  forecast  methodologies  is  to  create  a  set  of  closely -related,  physically 
realizable  initial  states  from  which  separate  forecast  integrations  can  begin.  As  the  forecast  in¬ 
tegrations  evolve,  the  ensemble  members  will  spread  (cluster)  about  uncertain  (relatively  certain) 
forecast  states.  The  forecast  ensemble  response  quantifies  the  uncertainty  in  each  forecast  state  in 
terms  of  the  spread  (i.e.  RMS  differences). 

The  initial  condition  concentration  of  MFS  ensemble  forecast  uncertainty  remains  in  the  mesoscale 
eddies  for  the  duration  of  the  MFS  re-forecast  experiments  run  to  date.  Figure  4  depicts  the  SSH 
spread  at  days  3  (Fig  4a)  and  10  (Fig  4b)  of  the  forecast  period.  The  largest  spread  occurs  in 
mesoscale  patterns  associated  with  the  Algerian  Current  in  the  Western  Basin,  and  with  a  bound¬ 
ary  current  off  North  Africa  in  the  Eastern  Basin.  The  Algerian  Current  instabilities  are  being 
analyzed  in  greater  detail  by  Bonazzi  et  al.  (2006). 

Animations  of  2  forecast  realizations  in  the  Gulf  of  Lyons  were  created  for  the  presentation 
at  ONR  HQ.  The  animations  demonstrate  significantly  different  timings  and  quantities  of  deep 
water  formation  in  the  pre-conditioned  region  of  the  Gulf.  This  suggests  a  sub-basin  to  general 
circulation  scale  implication  of  the  different  forcings  from  the  wind  realizations  drawn  from  the 
MFS-Wind-BHM  posterior  distribution.  The  sub-basin  scale  ocean  response  differences  are  due 
to  the  impact  of  different  mesoscale  evolutions,  driven  by  subtle  differences  in  the  wind.  Because 
only  one  ensemble  member  exhibits  earlier  and  more  abundant  deep  water  formation,  we  can 
characterize  this  particular  forecast  as  less  likely,  but  nonetheless  physically  consistent  within  the 
bounds  of  uncertainty  that  propagate  from  the  surface  wind  forcing. 

The  initial  condition  uncertainty  distribution,  and  the  retention  of  forecast  uncertainties  in  the 
mesoscale  eddies  is  a  significant  achievement  of  the  MFS  ensemble  forecast  system,  given  surface 
wind  forcing  realizations  from  MFS-Wind-BHM. 

MFS-Error-BHM  Plans  for  Year  2 

In  year  2  of  the  project  we  will  investigate  BHM  methods  for  modelling  the  evolution  of  error 
covariance  in  MFS.  In  the  MFS  ROOI  method,  the  full  covariance  matrix  P  is  replaced  by  a 
background  error  covariance  B.  The  forecast  background  error  covariance  matrix  is  decomposed 

Bf  =  d/1/2c-1D/1/2  (3) 
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where  is  a  diagonal  matrix  with  diagonal  elements  quantifying  the  forecast  error  variances,  and 
C  is  a  spatial  correlation  matrix. 

The  evolution  of  the  forecast  error  variances  is  ad  hoc,  we  let  =  NY)  '.  MFS-Error-BHM 
will  be  developed  to  model  N. 
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Figure  1:  Snapshot  of  surface  wind  vector  fields  for  the  Western  Mediterranean  from  MFS-Wind-BHM. 
Ten  realizations  from  the  posterior  distribution  are  shown  at  each  output  grid  point  (blue  vectors).  The 
ensemble  mean  vector  is  green  and  a  red  vector  denotes  the  ECMWF  analysis  vector  for  the  same  time.  For 
clarity,  wind  vectors  are  plotted  at  every  other  grid  point  in  the  MFS-Wind-BHM  prediction  grid. 
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Figure  2:  MFS  data  assimilation  and  forecast  cycle.  Full  ocean  data  assimilation  steps  are  run  every  14 
days  (blue  bars).  A  forecast  is  launched  at  the  end  of  each  assimilation  step  (red  arrow).  Forecasts  arc  also 
launched  on  days  between  data  assimilation  steps,  after  a  daily  simulation  run  (green  arrows),  wherein  only 
the  atmospheric  forcing  is  updated  (i.e.  ECMWF  analyses  replace  atmospheric  forecast  fields). 
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Figure  3:  MFS  ensemble  initial  condition  spread  (RMS  differences)  in  sea  surface  height  (SSH  in  cm;  top 
panel)  and  sea  surface  temperature  (SST  in  °C;  bottom  panel).  Initial  condition  differences  are  due  to  slight 
differences  in  wind  forcing,  for  the  data  14-day  data  assimilation  step,  from  across  the  10  realizations  of  the 
posterior  distribution  drawn  from  MFS-Wind-BHM.  In  both  SSF1  and  SST  the  differences  accumulate  in  the 
mesoscale  eddy  field. 
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Figure  4:  MFS  forecast  spread  in  SSH  (cm)  at  forecast  days  3  (top  panel)  and  10  (bottom  panel).  Ten  MFS 
forecasts  are  forced  by  different  realizations  from  the  MFS-Wind-BHM  posterior  distribution.  Forecasts 
are  launched  from  corresponding  initial  conditions  that  were  also  generated  by  BF1M  winds  in  the  data 
assimilation  step.  The  initial  condition  spread  that  was  concentrated  in  the  oceanic  mesoscale  (see  Fig.  3), 
persists  at  that  scale  throughout  the  10-day  forecasts. 
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